function [radVector, TwoCellVector] = cellPositionVector(cellNum,imfile,PrePost,dataflag)
% FUNCTION [RADVECTOR, TWOCELLVECTOR] = CELLPOSITIONVECTOR(CELLNUM, IMFILE...
% PREPOST, DATAFLAG) finds radial vector RADVECTOR, and direction vector
% TWOCELLVECTOR of CELLNUM. 
% written by Taro Kiritani, tarokiritani2008@u.northwestern.edu 12/14/2010

% find the radial vector.
radVector = findRadVector(cellNum, imfile, 1);

% find direction vectors.
TwoCellVector = getCellPosition(cellNum, PrePost);

% save radVector and TwoCellVector.
cellDistfile = ['C:\Data\Taro\CELLS\TK',cellNum,'\cellDist.mat'];
save(cellDistfile,'radVector','TwoCellVector')

end
    
function TwoCellVector = getCellPosition(cellNum, PrePost)
% FUNCTION TWOCELLVECTOR = GETCELLPOSITION(CELLNUM, PREPOST) returns
% direction vectors of cells. The vectors are calculated from either the .m
% file or picture files from the experiments.

cellPosition = readXYZfromMfile(cellNum);
% check if [x, y, z] location is stored in the experiment .m file.
vecSignal = [];
for posNum = 1:length(cellPosition)
    vecSignal(posNum) = isvector(cellPosition{posNum,1});
end

if sum(vecSignal) == 0 % coordinates not found in the .m file.
    TwoCellVector = getDirecVectorFromFigs(cellNum, PrePost)
else % coordinates found.
    for pairNum = 1:size(PrePost,1)
        TwoCellVector{PrePost(pairNum, 1),PrePost(pairNum, 2)} = cellPosition{PrePost(pairNum,2),1} - cellPosition{PrePost(pairNum,1),1};
    end
end
end
function radVector = findRadVector(cellNum, imfile, existFileFlag)
% FUNCTION RADVECTOR = FINDRADVECTOR(CELLNUM, IMFILE, EXISTFILEFLAG) either
% reads RADVECTOR from CELLDIST.MAT or manually measures radVector.
% written by Taro Kiritani, tarokiritani2008@u.northwestern.edu 1/17/2011

% read a mat file.
if existFileFlag == 1
    try
        load(['C:\Data\Taro\CELLS\TK',cellNum,'\cellDist.mat'])
    end
end

% manually measure.
if existFileFlag ~= 1 || not(exist('radVector','var'))
    % find the radial axis in 4x images.
    x4image = imread(['C:\Data\Taro\CELLS\TK',cellNum,'\images\',imfile]);
    figure;imagesc(x4image);colormap('gray')
    daspect([1 1 1]);set(gcf,'Position',(get(0,'ScreenSize')))

    % find the radial vector based on the 4x image loaded above. Note that
    % radVector is based on Euclidian coordinate system. The right and up
    % directions are positive in the 4x image. Also, the slice is assumed
    % to be flat and the radial axis is parallel to the focal plane.
    disp('first click pia, then white matter.');
    goSignal = 0;
    while goSignal ~= 1
        radVector = ginput(2);
        line(radVector(:,1),radVector(:,2));
        goSignal = listdlg('ListString',{'ok','select again'},'ListSize',[100 100]);
    end

    close
    radVector = radVector(1,:) - radVector(2,:);    
    radVector(2) = -radVector(2);
    radVector(3) = 0;
end

end
function cellPosition = readXYZfromMfile(cellNum)
% load .m file from the indicated experiment folder.
cd(['C:\Data\Taro\CELLS\TK',cellNum]);
mRead = cellstr(caseread(['TK',cellNum,'.m']));
cellPositionIndex = strfind(mRead, 'subregion');
cellPositionIndex = find(cellfun(@isempty,cellPositionIndex)==0);

for kk = 1:length(cellPositionIndex)
    eval(mRead{cellPositionIndex(kk),1});
end

cellNum

% read the positions of the cells from the .m file.
% cellPosition is based on Euclidian coorinate system. Note that y is 
% inverted in Rig 4.
cellPosition={'subregion1';'subregion2';'subregion3';'subregion4'};

for k = 1:length(cellPosition)
    try
        % there might be up to four cells. If some of the do not exist, 
        % corresponding elements in CELLPOSIRION is not over-written.        
        cellPosition{k,1} = eval(cellPosition{k,1});   
    end
    
    % convert string to xyz vector.
    [x yz] = strtok(cellPosition{k,1},{'(',',',')'});
    [y z] = strtok(yz, {'(',',',')'});
    z = strtok(z,{'(',',',')'});
    % check the format of cellPosition.
    [x, ' ' , y, ' ', z];
    [str2num(x),str2num(y),str2num(z)];
    cellPosition{k,1} = [str2num(x),-str2num(y),str2num(z)];         
end

end
function TwoCellVector = getDirecVectorFromFigs(cellNum, PrePost)

cd images
showImageList;
figSet = findobj(gcf,'Type','image');
for pairNum = 1:size(PrePost,1)
    % first check the direction vector exists or not. then check the opposite direction has been measured already. Ex.    
    % if cell1 -> cell2 has been measured, skip cell2->cell1   
    % measurement.
    try
        load(['C:\Data\Taro\CELLS\TK',cellNum,'\cellDist.mat'])
    end
    
    if not(exist('TwoCellVector', 'var'))
        TwoCellVector = cell(4);
    end
        
    if isempty(TwoCellVector{PrePost(pairNum, 1), PrePost(pairNum, 2)}) ...
            && isempty(TwoCellVector{PrePost(pairNum, 2), PrePost(pairNum, 1)})
        % neither the direction vector nor the opposite does not exist.
        
        % find rig #.
        cd(['C:\Data\Taro\CELLS\TK',cellNum]);    
        mRead = cellstr(caseread(['TK',cellNum,'.m']));    
        eval(mRead{3,1})    
        Rignum = str2num(Rig);

        % which cells to choose?
        disp(['Pre -> Post: ',num2str(PrePost(pairNum,1)),...
            ' to ',num2str(PrePost(pairNum,2))])
        
        % find scale factor.
        objLens = input('4x or 60x? (input 4 or 60): ');
        if Rignum == 2 && objLens == 4 
            scaleFactor = 2600;    
        elseif Rignum == 2 && objLens == 60
            scaleFactor = 177;        
        elseif Rignum == 4 && objLens == 4
            scaleFactor = 2901;
        else
            scaleFactor = 204;    
        end

        % show a picture.
        selectedImg = input('which fig?: ');    
        figure;    
        selectedImg = get(figSet(length(figSet) - selectedImg + 1),'CData');    
        imagesc(selectedImg)    
        colormap('gray')    
        pause    
        set(gcf,'Position',(get(0,'ScreenSize')))

        % get TwoCellVector
        DirecVec = ginput(2) * scaleFactor / size(selectedImg,2);
        DirecVec = diff(DirecVec);
        DirecVec(2) = -DirecVec(2);
        DirecVec(3) = 0;
        TwoCellVector{PrePost(pairNum, 1), PrePost(pairNum, 2)} = DirecVec;
        close;
    elseif isempty(TwoCellVector{PrePost(pairNum, 1), PrePost(pairNum, 2)})
        % the opposite vector exists.
        TwoCellVector{PrePost(pairNum, 1), PrePost(pairNum, 2)} = ...
            -TwoCellVector{PrePost(pairNum, 2), PrePost(pairNum, 1)};
    end
end 
        % the long edge of field of view is 204 um in rig 4.
         % the long edge of field of view is 2901 um in rig 4.
         % the long edge of field of view is 177 um in rig 2.
         % the long edge of field of view is 2600 um in rig 2.
end